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Abstract 

Context. In the past, the isolated, radio-quiet neutron star RX J0720.4— 3125 showed variations in its spectral parameters 
(apparent radius, temperature of the emitting area and equivalent width of the absorption feature) seen in the X-ray 
spectra, not only during the spin period of 8.39 s, but also over time scales of years. New X-ray observations of 
RX J0720.4— 3125 with XMM-Newton extend the coverage to about 7.5 years with the latest pointing performed in 
November 2007. Out of a total of fourteen available EPIC-pn datasets, eleven have been obtained with an identical 
instrumental setup (full frame read-out mode with a thin filter), and are best suited for a comparative investigation of 
the spectral and timing properties of this enigmatic X-ray pulsar. 

Aims. We analysed the new XMM-Newton observations together with archival data in order to follow the spectral and 
temporal evolution of RX J0720.4-3125. 

Methods. All XMM-Newton data were reduced with the standard XMM-SAS software package. A systematic and 
consistent data reduction of all these observations is warranted in order to reduce systematic errors as far as possible. 
Results. We investigate the phase residuals derived from data from different energy bands using different timing solutions 
for the spin period evolution and confirm the phase lag between hard and soft photons. The phase shift in the X-ray 
pulses between hard and soft photons varies with time and changes sign around MJD=52800 days, regardless of the 
chosen timing solution. The phase residuals show a marked dependence on the energy band, and possibly follow a cyclic 
pattern with a period of ~9-12 yrs. We find that an abs(sine) dependence provides a better fit to the residuals with 
respect to a simple sine wave, although in both cases the reduced is not completely satisfactory. We compared the 
model fit to phase residuals derived from EPIC-MOS and Chandra (HRC and ACIS) data restricted to the hard energy 
band (400-1000 eV) to take into account the different energy responses of these instruments. 

Conclusions. The new data are not in contradiction with RX J0720.4— 3125 being a precessing neutron star but do not 
provide overwhelming evidence for this picture either. 
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1. Introduction 

The radio quiet, isolated neutron star RX J0720.4— 3125 
was discovered in ROSAT all-sky survey data 
(jHaberl et aL 1997| . It belongs to a group of seven 
neutron stars with similar properties often called the 
magnificent seven, hereafter M7. The X-ray spectra are 
well represented by a blackbody model (with absorption 
features in some cases) without any additional non-thermal 
component. For the two brightest stars (RX J1856. 5-3754 
and RX J0720.4— 3125) the trigonometric parallax was 
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measured (Kaplan et al., I2007P which enables us to 
estimate the size of the emission region. Six of them are 
X-ray pulsars with relatively long periods between 3.5 s 
and 11.4 s (for a review see, e.g., Haberl, I2007p . From 
the spin period, period derivative (P) and the energy 
of absorption lines one can independently estimate the 
magnetic field strength on the surface of the neutron star. 
Both methods consistently suggest strong magnetic fields 
in excess of ~ 10^"^ Gauss (Haberl, 2005) . The properties of 
the M7 make them of special interest for the investigation 
of cooling models of neutron stars (see Page & Reddy, 
12006| Page et al. l2006l and Blaschke & Grigorian, [2007|l 
and for constraining of the equation of state of matter at 
nuclear densities (Triimper et al.. i2004|) . 
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RX J0720.4— 3125 shows a long term variation in its 
spectral parameters (see de Vries et al., '2004': "Vink et al., 
[MHl and Haberl et al., [30(36: hereafter H06), i.e. its tem- 
perature, size of the emitting area and equivalent width of 
the absorption feature. The variations derived from data 
covering 5.5 years were consistent with a period of about 
(7.1 ± 0.5) yr, when assuming a sinusoidal shape of the 
variations (H06). Kaplan & van Kerkwijk (2005J (here- 
after KvK05) performed a coherent timing analysis for 
RXJ0720.4-3125 by including all data (ROSAT, XMM- 
Newton and Chandra) available at that time. This enabled 
them to determine the period more accurately then pre- 
vious studies (Cropper et al., 120011 Zane et al.. 12002] and 
Kaplan et al.. [3Ug3|l . In their Table 2 KvK05 provide two 
different values for a constant P: 7.016 x 10""'^''' s/s for 
Chandra data only and 6.983 x lO^^'* s/s using all data. 
By applying these timing solutions to derive phase binned 
light curves, KvK05 obtained phase residuals. These phase 
residuals, as shown by H06, can be modeled by a sinu- 
soidal function with a period of (7.7 ± 0.6) yr. This is 
consistent with the (7.1 ± 0.5) yr period of the spectral 
variations. As a possible explanation for the behaviour of 
RX J0720.4-3125, de Vries et al. (|20M) and H06 suggested 
precession of the neutron star, van Kerkwijk et al. (20(37| 
hereafter vK07) derived a new P = 6.957 x 10"^** s/s in- 
cluding new XMM-Newton and Chandra data, but still ob- 
tained sinusoidal-like phase residuals. KvK05 and vK07 de- 
rived x^=l for their timing solutions only after adding an 
artificial uncertainty to the time of arrivals (see KvK05 
and vK07 for more details). Furthermore, they add a glitch 
model to the spin-down term and claim the occurrence of 
the glitch around MJD=52800 days to explain the timing 
and spectral changes of this neutron star. Equivalent to the 
glitch model vK07 found a period w4.3 yr but that is not 
consistent with the period «7.5 yr found by H06. On the 
other hand, the data set used in H06 (EPIC-pn) is quite 
different to that used in vK07 (also including EPIC-MOSl 
& M0S2 and Chandra HRC-S/LETG). 

This paper includes new data from the recent XMM- 
Newton (Jansen et al. . 1200111 observations which are part of 
an ongoing monitoring programme, and Chandra archival 
data. We discuss our results as new evidence for a cyclic 
behaviour of RX J0720.4— 3125 based on timing solutions 
from KvK05 and vK07. 



2. X-ray observations and timing analysis 

Since the work of Haberl (I2007P was published, four new 
XMM-Newton observations have been performed between 
May 2006 and November 2007 at six month intervals (satel- 
lite revolutions 118 1 126 5, 1356 and 1454). All new EPIC- 
pn (Striider et al. I200ip observations were executed with 
the same instrument setup in full frame mode and with 
thin filter (see Table [l] for a summery of the observations) . 
Inclusion of the latest pointings brings the number of avail- 
able EPIC-pn datasets to fourteen. These, in conjunction 
with the MOS (Turner et al.. lM?T|) data (see Table [l] for 
details on the instrumental setup), have been used for the 
source timing analysis. All observations were reduced with 
the XMM-Newton Science Analysis System (SAS) version 
7.1.0. 

To derive the phase residuals from a certain timing so- 
lution, we binned the light curve in 40, 100 and 400 phase 



bins. The number of phase bins did not influence the result. 
Considering the time resolution of the full frame (2.6 s) and 
small window mode (0.3 s) we divide the MOS light curves 
in 10 and 40 phase bins, respectively. Phase residuals from 
all observations are obtained fitting a sinusoidal function to 
the pulse profile. The uncertainties of the phase residuals 
are derived from the uncertainties of the sinusoidal fit. All 
errors on the phase residuals in this paper correspond to 
Icr. 

In the analysis of KvK05 the phase residuals for the 
last observations show an increasing trend while in the 
case of vK07 (with more data included) the trend is re- 
versed. Therefore, the timing solution for RX J0720.4— 3125 
strongly depends on the behaviour of the source in the fu- 
ture and all solutions obtained so far must be regarded 
as preliminary. Hence, in the following we examine phase 
residuals for two timing solutions provided by KvK05 and 
vK07 (in both cases their "all data solution"). 

To reduce systematic effects due to the changing low- 
energy response of the EPIC-MO S cameras we use only 
events from the 400-1000 eV band for MOS (the 120-1000 
eV band is used for pn). The phase residuals using the co- 
herent timing solution of KvK05 are shown in Fig. [l] As 
already mentioned in Haberl ()2007p . systematic differences 
between phase residuals from different instruments can be 
seen, which are probably caused by the energy dependent 
pulse profile. To verify this, we investigated the phase resid- 
uals for two different bands (soft: 120-400 e'V and hard: 400- 
1000 e'V) using EPIC-pn data, and the results are shown 
in Fig. |2] The phase residuals exhibit a clear energy depen- 
dence while the phase shift between the two bands varies 
with time and changes sign around MJD=52800 days. As 
an example of a relatively large phase shift between hard 
and soft bands we show the pulse profile for satellite orbit 
986 in Fig.[3j The phase shift is clearly visible and the trend 
is much the same as in Fig. [2] if we add or subtract the error 
of P provided in KvK05. 

'We also applied the timing solution from KvK05 to 
archival Chandra HRC-S/LETG (zeroth order) and ACIS 
data (for details of the observations see KvK05 and vK07) . 
For observations 368, 369, 745, 2771, 2772, 4667, 4668, 4669, 
4670, 4671, 4672, 7177, 7243 and 7245 we use 10 phase 
bins and 40 phase bins for the other observations because 
of the different time resolution and number of photons. 
Considering the energy dependence of the phase residuals 
we restrict the energy band for the analysis of the ACIS 
data to the hard band (400-1000 e'V). 'We show all phase 
residuals from Chandra and XMM-Newton observations in 
Fig. |4] Comparison with Fig. [l] shows that the agreement 
between the different instruments is much improved. 

As in the EPIC-pn observations before (see Haberl, 
l2007i H06 and de 'Vries et al., l2004ll . the pulse profile can 
be approximated by a sinusoid. In Fig. [5] we show the pulse 
profiles for the last four EPIC— pn data sets together with 
the hardness ratios. The hardness ratios HR = H/S are 
anti-correlated with intensity. The pulse phases were de- 
rived by applying the "all data" solution of KvK05 (2005). 
The light curves look the same, except for absolute phase of 
course, using the timing solution of vK07. As noted before, 
the profiles are not exactly sinusoidal and look somewhat 
skewed. This is particularly true for the latest observations. 
To explore how much the behaviour of the residuals is in- 
fluenced by the choice of the pulse shape we also tried a 
sawtooth proflle, in which the location and the height of 
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Figure 1: Phase residuals of RX J0720.4— 3125 derived from 
the timing solution from Kaplan & van Kerkwijk (|2005p . 
Red dots: refer to EPIC-pn data, blue dots: to EPIC-MOSl 
& M0S2 observations, black dots: reproduced from the 
original solution of Kaplan & van Kerkwijk (|2005p . 

the minimum are free to vary. Results for the phase residu- 
als in the two bands (see Fig. [6]) are similar to those derived 
for a sinusoidal profile. This supports the fact that energy- 
dependent phase residuals are indeed present and do not 
depend on the chosen shape of the pulse profile. We note 
that the jd.o.j values of the sine fits shown in Fig. are 
smaller than those of the sawtooth fits {x^ /d.o.f = 2.94, 
4.06, 3.70, 2.50 in the order of increasing revolution num- 
ber) so that a sine wave provides a better interpretation of 
the data, although the quality of the fits is still not satis- 
factory. For all fits d.o.f — 37. 

The timing solution of vK07 with a glitch at 
MJD=52817 days produces phase residuals three times 
larger than those seen in Fig. [2] with a large scatter and 
will not be considered any further. The pulse profiles for 
EPIC— pn satellite orbit 622 and 711 in both energy bands 
in Fig. [7] are an example of the large phase shifts from one 
observation to the next by applying this solution. 

3. Modeling the long term evolution of the timing 
properties 

From the first XMM-Newton observations of 
RXJ0720.4-3125 Cropper et al. (PUUT|| reported a 
phase shift between X-ray intensity and hardness ratio 
(HR, hard and soft band count ratio). Furthermore, de 
Vries et al. p004p found that the phase shift changed 
sign between the year 2000 (the HR leads the intensity 
profile) and 2003 (the HR lags behind the intensity). As 
can be seen in Fig. |2] this occurred around MJD=52800 
days. Since then, the hard band emission lags behind the 
soft band emission and the shift seems to decrease again 
from a maximum of ~ 0.1. The phase lag dependence 
on the energy bands clearly shows that using data from 
instruments with different energy-dependent sensitivity 
may introduce systematic errors in the timing analysis. 
The EPIC-pn instrument has a high low-energy response 
relative to the hard band, therefore Chandra ACIS and 







I II I I I 




1 1 


T T 


- I ' 




I " ■ 1 1 




\ I 

] I 


I 















_0 2l ^ ^ ^ ^ ^ 1 

5.15 5.2 5.25 5.3 5.35 5.4 5.45 



MJD [days] 



Figure 2: Phase residuals of RX J0720.4— 3125 derived from 
EPIC-pn data with the timing solution from Kaplan & van 
Kerkwijk PK?51) in the soft band (120 - 400 eV, red) and 
hard band (400 - 1000 eV, black). 
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Figure 3: Phase shift of hard (black) and soft (red) band 
of RX J0720.4-3125 derived from EPIC-pn data (sateUite 
orbit 986) with the timing solution from Kaplan & van 
Kerkwijk ()2005p . Errors are Poissonian. The y-axis denotes 
the relative counts per phase bin. 



EPIC-MOS are expected to be closer to the EPIC-pn 
results restricted to the hard band. 

As in previous work, we start modeling the phase resid- 
uals using a sine wave function (Fig.|8] top). To avoid sys- 
tematic effects we use only the EPIC-pn data (120-1000 eV) 
for our timing analysis. The different period derivatives, P, 
given in KvK05 and vK07 depend on which and how many 
data points and from which instruments they were used. 
We already mentioned the different energy response of the 
different instruments. P depends on how much of the cy- 
cle, shown in Fig. [l] is covered with observations and can 
be determined more accurately only in the future, when 
at least one full cycle of the variation of the phase resid- 
uals is observed, assuming a 10-12 yr period as discussed 
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Table 1: XMM-Newton EPIC observations of RX J0720.4-3125. 



Orbit 


Observation 


Date 




MJD 


EPIC instrument setup^ 


Duration [s] 










[days] 


pn 


MOSl 


M0S2 


(EPIC-pn) 


78 


0124100101 


2000 


May 13 


51677 


FF thin 


FF thin 


SW thin 


62498 


175 


0132520301 


2000 


Nov. 21-22 


51869 


FF medium 


SW thin 


- 


26118 


533 


0156960201 


2002 


Nov. 6-7 


52584 


FF thin 


FF thin 


FF thin 


28373 


534 


0156960401 


2002 


Nov. 8-9 


52586 


FF thin 


FF thin 


FF thin 


30173 


622 


0158360201 


2003 


May 2-3 


52761 


SW thick 


FF thin 


FF thin 


72796 


711 


0161960201 


2003 


Oct. 27 


52939 


SW medium 


FF thin 


FF thin 


43013 


815 


0164560501 


2004 


May 22-23 


53147 


FF thin 


FF thin 


FF thin 


43951 


986 


0300520201 


2005 


April 28 


53488 


FF thin 


SW thin 


SW thin 


51435 


1060 


0300520301 


2005 


Sep. 23 


53635 


FF thin 


SW thin 


SW thin 


51135 


1086 


0311590101 


2005 


Nov. 12-13 


53686 


FF thin 


SW thin 


SW thin 


37835 


1181 


0400140301 


2006 


May 22 


53877 


FF thin 


SW thin 


SW thin 


20035 


1265 


0400140401 


2006 


Nov. 05 


54044 


FF thin 


SW thin 


SW thin 


20035 


1356 


0502710201 


2007 


May 05 


54225 


FF thin 


FF thin 


FF thin 


20035 


1454 


0502710301 


2007 


Nov. 17 


54421 


FF thin 


FF thin 


FF thin 


23059 



'^'Read - out mode and filter; FF; Full Frame; SW; Small Window. 




Figure 4: Phase residuals of RX J0720. 4-3125 derived from 
the timing solution from Kaplan & van Kerkwijk (j2005p . 
Black dots: EPIC-pn, red dots: EPIC-MOSl small window 
and thin filter, red stars: EPIC-M0S2 small window and 
thin filter, blue dots: EPIC-MOSl full frame and thin filter, 
blue stars: EPIC-M0S2 full frame and thin filter (all in 
hard band); open circles show phase residuals from Chandra 
HRC (no energy selection) and ACIS (hard band). 



below. Therefore, we introduced an additional linear term 
y=mt-|-n, as can be seen in Fig. [8] for the timing solution of 
KvK05 (middle) and that of vK07 with constant spin-down 
without glitches (bottom). 

The two maxima in the phase residuals (Fig. [8| do not 
necessarily have equal heights. By adding the linear term, 
we automatically take into account the different values of 
P and a possible different amplitude of the two maxima 
covered by the observations. For the three cases described 
above, we derive periods of (5.58 ± 0.33) yr, (7.20 ± 1.50) yr 
and (5.59 ± 0.57) yr, respectively (here and below, errors 
on the periods correspond to la). The values are consis- 
tent within their errors and in particular the first and last 
value are nearly identical. This is expected because a pos- 
sible period should not significantly depend on the actual 



timing solution (as long as this is correct to first order). It 
should be noted, however, that the fits with a sine wave are 
formally not acceptable, given the derived values. The 
unacceptable fit can be partly explained by the two data 
points from MJD=52584 and 52586 days. It suggests that 
the sine wave is not a good representation of the data. We 
excluded the data point from MJD=52584 days and de- 
rived periods of (6.54 ± 0.14), (6.55 ± 0.46) and (5.73 ± 
0.43) yr, respectively, for the three cases, while the values 
for xVd-o.f are 17.16/9, 16.90/8 and 35.12/8, respectively. 
Only the fit with the data shown in Fig. [8](top), excluding 
the data point at MJD=52584 days, is formally acceptable. 
The reason for the outlier at this epoch is unclear. 

In the first observation the soft photons lag behind the 
hard photons, while in the observations after MJD=52800 
days the opposite is true (see Fig. |2|) . The phase lag in the 
first observation is supported by Fig. 2 of de Vries et al. 
(|2004p . which is independent of a model fit to the pulse pro- 
file. The energydependent heights of the two phase residual 
maxima (Fig. 2]) show that the actual period of the mod- 
ulation could be twice the sine period. The energy depen- 
dence suggests a model of emission originating from two hot 
spots with different temperature. In the case of a precession 
model, for one half of the precession period the soft emission 
(cooler pole) would precede the hard emission and for the 
other half it would be the other way around. The cooler cap 
would be mostly visible for half the precession cycle, while 
in the other half the hotter cap would dominate emission. 
If the emission arises in different poles this suggests that 
the profile should be that of an abs(sine), similar to that 
used by vK07 and van Kerkwijk & Kaplan (I2007P but with 
a different interpretation. They found a period R:i4.3 yr as a 
solution equivalent to the glitch model. Unfortunately van 
Kerkwijk & Kaplan (2007) do not provide errors, and we 
cannot check if our value of the period ^5.5 yr is consistent 
with theirs. The difference of the values may be explained 
by the four new EPIC-pn observations we add. 

We interpret the turn-around in the phase residual evo- 
lution (near MJD=52800 days) not as a glitch, but as 
"switching" of the emission to the predominantly visible 
pole. As in the fit with the normal sine function, we again 
add the linear term to the abs(sine). To reduce systematic 
effects due to the use of different filters we use the EPIC- 
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Figure 5: Folded EPIC-pn pulse profile (120-1000 eV, 40 phase bins, d.o.f = 37 in all cases) with sinusoidal fit and hardness 
ratio HR=H/S (bands described in the text) of the last four observations. Pulse phases are obtained by applying the "all 
data" solution of Kaplan & van Kerkwijk (|2005p . Errors are Poissonian. I/Imax denotes the relative counts per phase 
bin. 
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pn hard band (400-1000 eV) data. We use the hard band 
in order to be able to include data from other instruments, 
although statistics in the soft band would be preferable. 
For the fit we further increase the size of the error of the 
phase residual in revolution 622 by a factor of two (model 
B1-B4). This observation was performed with the thick fil- 
ter, which has a significantly lower transmission in the soft 
energy band compared to the medium and thin filter. This 
causes a systematic shift of this phase residual to larger val- 
ues and would introduce systematic errors in the fit results. 
Excluding this observation completely (model A1-A4), we 
obtained a period which is 2 - 4.5 years longer (see Tablel2]). 
To get an estimate of the systematic error of the period, 
which is likely larger than the statistical error derived from 
the fit, we also applied a doubled error for all observations 
which were not performed with the thin filter, i.e. revolu- 
tions 175, 622 and 711 (model C1-C3). The results of the 
fitting procedures are listed in Table [2] and plotted in Fig. 
[9) Except for model A, the advantage of introducing an ad- 
ditional slope is clearly visible in the x^/d.o.f. as compared 
to the poor fits without a slope. 

Applying a sine instead of the abs(sine) models Bl, B2 
and B3, we derive x^/d.o.f.=3.31, 3.67 and 4.85, compared 
to 12.25 (top), 2.43 (middle) and 4.99 (bottom) in Fig. |9] 
The fit of model Bl is much better with a sine instead of an 
abs(sine) but a better approximation would be an abs(sine) 
with an additional slope, which gives the best x^=2.43 (B2) 
in this case for the timing solution of KvK05. The two 
values for the timing solution of vK07 are not significantly 
different. 



For a consistency check we show in Fig. 10 the abs(sine) 
model fit shown in Fig. |9] (middle) together with the phase 
residuals obtained from EPIC-pn, EPIC-MOS and ACIS 
(all hard band) and HRC observations. We excluded in this 
Figure Chandra data points with errors > 0.03. 

4. The long term spectral evolution 

To avoid systematic errors due to cross-calibration uncer- 
tainties between the different instruments we restrict our 
spectral analysis to the EPIC-pn data. Moreover we use 
only the eleven EPIC-pn observations of RX J0720.4-3125, 
which were done in full frame mode with the thin filter. 
We extracted the pulse-phase averaged spectra from cir- 
cular source and background regions with 1' in diameter, 
selecting single-pixel events (PATTERN = 0). The spec- 
tra were binned to obtain at least 100 photons in each bin. 
The fits were performed with XSPEC12 and restricted to 
energies between 0.16 keV and 1.5 keV. We used the black- 
body model plus an additive Gaussian line, which repre- 
sents the absorption feature (negative normalisation). The 
energy resolution of the instrument is not sufficient to de- 
termine the shape of the absorption feature. A Gaussian 
profile applied multiplicative ly yields equally good fits. The 
feature might also consist of several unresolved lines or have 
the shape of an absorption edge. Here and below, the errors 
on the spectral parameters correspond to the 2.7ct confi- 
dence level. 

We linked all model parameters for the spectra of revo- 
lution 533 and 534, because the two observations were only 
two days apart. The line energy, the line width a and the 
hydrogen column density Njj were linked for all spectra, 
i.e. we assume that they do not vary with time. All spec- 
tra were fitted simultaneously. The best fit parameters are 



(301 ± 3) eV for the Hue energy, cr = (77 ± 2) eV and 
Nh = (1.04 ± 0.02)-1020 cm-2. For all spectra together we 
obtain x^/d.o.f=1.29 with 1767 degrees of freedom. The 
blackbody temperature kT, the blackbody normalisation 
(used to calculate the size of the emitting area, assuming 
a distance of 300 pc) and the line equivalent width EW 
were free parameters and their best fit values are listed in 
Tab le [3| The time evolution of the parameters is drawn in 



Fig. [n[ We also included kT and EW of the three EPIC- 
pn observations done with different observing mode and 
different instrument setup in this figure, which were not 
used for this fit. Because of the cross-calibration uncertain- 
ties mentioned above, the blackbody radii obtained from 
these observations are not reliable and thus not included in 
Fig. 
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bottom. 

As can be seen in Fig. |11| the sine wave form is not a 
good approximation for the new data. We only can conclude 
that, if the variation in the spectral parameters is cyclic, 
the period is larger than the time span of the EPIC-pn 
observations. Currently it is not possible to predict whether 
it is periodic or not and if it is periodic, what the shape of 
the variation might be. 

Following H06 we computed pulse phases for each de- 
tected event using the timing solution of KvK05. We di- 
vided one pulse period into five phase intervals and ex- 
tracted spectra from each phase interval for each observa- 
tion. As in the analysis of the phase-averaged spectra, we 
performed a joint fit to the 55 spectra with the same param- 
eters linked (x^/d.o.f=1.16 with 6452 degrees of freedom). 
The change of kT and EW during the pulse period for each 
of the observations is shown in Fig. [T2j As seen before, 
the parameters evolve counter-clockwise and the long term 
trend of the parameters, as seen in Fig. 
pulse phases. 



5. Discussion 

We present new XMM- Newton observations of the isolated 
neutron star RX J0720.4— 3125 performed during our half- 
year monitoring program in 2006 and 2007. Combining the 
new data with the previous XMM-Newton data presented 
by H06 and with Chandra data (KvK05, vK07), we fol- 
lowed the spectral and timing behaviour of the source. As 
shown by KvK05 and vK07, a phase-coherent timing analy- 
sis assuming a constant spin-down of the pulsar yields rela- 
tively large phase residuals. Different models have been fit- 
ted to these phase residuals invoking a more sudden change 
in properties (vK07) on one hand or a periodic behaviour 
(H06 and vK07) on the other hand. Physically, these effects 
could be interpreted, e.g., by a "glitch" or by precession of 
the neutron star, respectively. 

The analysis of all available data from the different in- 
struments on XMM-Newton and Chandra yields system- 
atic shifts between the phase residuals (Fig.[T]). As already 
noted by Haberl p007p this is most likely caused by intrin- 
sicaly energy dependent phase residuals and this must be 
accounted for when using different instruments of different 
energy response. This is confirmed by our timing analysis 
using events from the EPIC-pn data in the soft (120 - 400 
eV) and hard (400 - 1000 eV) energy band separately, which 
yield systematic phase shifts (Fig. ^1. Moreover, the phase 
shifts between soft and hard photons vary with time. In 
order to still combine the results from the different instru- 
ments, we restricted our analysis to the hard energy band 
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Table 2: Periods derived from the fits of the phase residuals with a model consisting of abs(sinc)+mt+n for different 
timing solutions. The plots for models Bl-3 are shown in Fig. [91 All errors correspond to Icr. 



Model 


Timing 


period 


m amplitude 


offset (n) 




d.o.f 


xVd.o.f 




solution 


[yr] 


[phase res./s] [phase res.] 


[phase res.] 








revolution 622 excluded 


Al 


KvK05 


15.21 ± 0.65 


0.163 ± 0.068 


-0.079 ± 0.017 


24.51 


9 


2.72 


A2 


KvK05 


15.210 ± 0.011 


7.45 -10"" 0.15 ± 0.36 


-0.41 ± 3.40 


22.75 


8 


2.84 


A3 


vK07 


12.488 ± 0.055 


-1.10 -10"" 0.125 ± 0.034 


-0.01 ± 0.41 


64.49 


8 


8.06 


A4 


vK07 


10.693 ± 0.059 


0.065 ± 0.028 


-0.015 ± 0.015 


61.53 


9 


6.84 


doubled error of revolution 622 


Bl 


KvK05 


10.711 ± 0.058 


0.097 ± 0.050 


-0.028 ± 0.035 


122.45 


10 


12.25 


B2 


KvK05 


11.504 ± 0.072 


39.21 -10"" 0.099 ± 0.029 


-1.83 ± 0.40 


21.87 


9 


2.43 


B3 


vK07 


10.546 ± 0.069 


3.44 -lO"" 0.095 ± 0.034 


-0.20 ± 0.42 


44.89 


9 


4.99 


B4 


vK07 


9.656 ± 0.023 


0.057 ± 0.026 


-0.012 ± 0.017 


51.38 


10 


5.14 


doubled error of revolution 175, 622 and 711 


CI 


KvK05 


9.07 ± 0.14 


0.058 ± 0.068 


0.003 ± 0.042 


138.33 


10 


13.83 


C2 


KvK05 


9.376 ± 0.037 


45.17-10-" 0.087 ±0.035 


-2.09 ± 0.48 


18.45 


9 


2.05 


C3 


vK07 


9.009 ± 0.037 


2.48 -10"" 0.090 ± 0.040 


-0.15 ± 0.53 


32.61 


9 


3.62 


C4 


vK07 


9.017 ± 0.038 


0.076 ± 0.031 


-0.025 ± 0.019 


45.00 


10 


4.50 



0.15 



0.05- I I i I \ ^ - 

i - I * , ' I I- 

CD .. J- 

-0.05- 4 ? 

-1= T 
Q. " I 

-0.1- \ 

-0.15- J 

-0.2 ' ' ' 1 ' 
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Figure 6: Phase residuals of RX J0720.4— 3125 derived from 
the timing solution from Kaplan & van Kerkwijk (2005), 
similar to those in Fig. [2] but fitting a sawtooth instead of 
a sinusoidal to the pulse profile. 



for the three EPIC instruments and for ACIS. The negli- 
gible spectral capabilities of the HRC prevented an energy 
selection and we excluded the ROSAT data completely due 
to the large errors and alias periods caused by the inter- 
rupted observations. This brought the phase residuals from 
the various instruments in much closer agreement (Fig. |4]) . 

Because of the systematic shifts of the pulse phases on 
time scales of years, which become visible as phase resid- 
uals with respect to the model with constant spin-down 
of the neutron star, the solutions for P strongly depend 
on the time interval used for the analysis. Only when the 
observations eventually cover a period significantly longer 
than the time scale of the phase variations - which is cur- 
rently not the case - a reliable value for P can be derived. 
This is true if the phase variations are periodic, but also 
for a sudden event with a long (>10 yrs) relaxation time. 
To overcome this problem, we included a linear term into 
our long term modeling of the phase residuals and inves- 



Table 3: All EPIC-pn observations in full frame mode with 
the thin filter of RX J0720.4-3125 with the derived varia- 
tions of the pulse averaged spectra. Observations done in 
different observing modes, which were not used for the fit 
are in italic. The radius obtained from observation 711 is 
not used at all (see explanation at the end of section 5). All 
errors correspond to 2.7 cr confidence level. 



Orbit 


MJD 


kT 




EW 




radius 






[eV] 




[eV] 




[km] 


78 


51677 


86.5 


± 


0.4 


-4.0 ± 4.0 


4.79 ± 


0.09 


175 


51870 


86.0 


± 


0.7 


+5.6 ± 


70 


4.90 ± 


0.14 


533/534 


52585/7 


88.4 


± 


0.4 


-16.8 ± 


3.3 


4.70 ± 


0.07 


711 


52940 


92.0 


± 


0.6 


-64.7 ± 


5.4 






815 


53148 


94.6 


± 


0.5 


-58.6 ± 


4.2 


4.29 ± 


0.08 


986 


53489 


94.3 


± 


0.4 


-57.8 ± 


3.8 


4.35 ± 


0.07 


1060 


53636 


93.7 


± 


0.5 


-55.0 ± 


3.8 


4.37 ± 


0.08 


1086 


53687 


93.1 


± 


0.5 


-55.9 ± 


3.8 


4.50 ± 


0.08 


1181 


53877 


93.1 


± 


0.6 


-51.2 ± 


4.8 


4.44 ± 


0.10 


1265 


54045 


92.7 


± 


0.6 


-53.6 ± 


4.7 


4.53 ± 


0.11 


1356 


54226 


92.4 


± 


0.6 


-44.1 ± 


4.9 


4.45 ± 


0.11 


1454 


54421 


91.8 


± 


0.6 


-44.6 ± 


4.6 


4.50 ± 


0.10 



tigate two solutions for P as inferred by KvK05 (all data 
solution) and vK07 (all data solutions without glitches). 
However, we emphasize that both solutions should be re- 
garded as first approximations and that any future timing 
analysis combining data from different instruments needs 
to be restricted to narrower energy bands. 

Using an abs(sine) function to model the variations of 
the phase residuals, we obtain periods of 9-15 years with 
the most reliable values around 9.3-11.5 years for the two 
fits with the lowest (Table p|. The observations still 
cover less than two "humps" of the abs(sine) function and 
more data are needed to ultimately verify that the timing 
behaviour is indeed cyclic and to determine the period more 
precisely. The energy dependence of the phase shifts (Fig. |2| 
clearly demonstrates that the long term period covers two 
humps of the abs(sine). This was not taken into account 
previously and is at least partly the reason for the shorter 
periods determined in the past (H06, vK07). 

For the spectral analysis of RX J0720.4— 3125 it is even 
more important to avoid systematic differences due to cross- 
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Figure?: Pulse profiles of satellite orbit 622 (black) and 711 
(red) in soft (upper panel) and hard (lower panel) bands as 
an example for the large phase shift applying the timing 
solution from vK07 with a glitch. The phase residuals are 
much larger than for the other timing solutions (discussed 
in the text) and the pulse profile looks more skewed. Error 
bars denote Poissonian errors. On the vertical axis are the 
counts are normalized to their maximum value. 



calibration uncertainties between the various instruments. 
Therefore, we restricted our spectral analysis to the EPIC- 
pn data taken in full frame CCD readout mode and with 
a thin optical blocking filter. As in HOG we use a simple 
absorbed blackbody model including an absorption feature 
which we assume to be of Gaussian shape. The new XMM- 
Newton observations show that the inferred spectral param- 
eters continue to change, but do not follow the sinusoidal 
evolution suggested by the previous subset of data used in 
HOG. Temperature, kT, the apparent radius of the emitting 
area and the equivalent width, EW, of the absorption line 
changed more slowly over the years 200G and 2007 than an 
extrapolation of the sine model would predict. The spec- 
tral data do not allow us to draw any conclusions about 
periodicity. 

The phase lags between soft and hard emission (Fig. [2| 
suggest the existence of at least two hot regions on the neu- 
tron star surface with somewhat different emission charac- 
teristics. These hot spots might be associated with the polar 
caps of the neutron star in a magnetic field. Two hot spots 
have also been suggested to explain the rotational varia- 
tions seen in the light curves of RBS 1223 - another M7 
star - in different energy bands ( |Schwope et al. 2005D . To 
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Figure 8: Phase residuals of RX J0720.4-3125 derived from 
the EPIC-pn data using the timing solution from Kaplan 
& van Kerkwijk (2005) fitted with a pure sinusoidal model 
(top) and a sinusoidal model including a linear term (mid- 
dle). The fit with the latter model for the timing solu- 
tion of van Kerkwijk et al. (|2007p is shown in the bottom 
panel. The periods are (5.48 ± 0.33), (7.20 ± 1.50) and 
(5.59 ± 0.57) yr, respectively, while the values for x^/d.o.f 
are 59.74/10 (top), 43.05/9 (middle) and 57.25/9 (bottom 
panel) . 
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Figure 9: Phase residuals of RX J0720.4— 3125 derived from 
the EPIC-pn 0.4-1.0 keV data. The error of the phase resid- 
ual from revolution 622 (around MJD=52800 days) was in- 
creased by a factor of two (see text) . Shown are the best fit 
models based on different abs(sine) functions and different 
timing solutions as summarized in Table |2] The derived pe- 
riods are (10.711 ± 0.058) yr (model Bl, top), (11.504 ± 
0.072) yr (B2, middle) and (10.546 ± 0.069) yr (B3, bot- 
tom). 



Figure 10: The phase residuals (only hard band) shown in 
Fig. [4] together with model B2 derived from the EPIC-pn 
(hard) data shown in Fig. [9] (middle). Phase residuals from 
ROSAT and Chandra with errors larger than 0.03 (phase 
residuals) are excluded. 



explain the different energy dependence of the two peaks are 
seen in the pulse profiles of RBS 1223, these authors suggest 
two hot spots with different size and temperature which are 
not anti-podal. A similar geometry can in principle also ac- 
count for the spectral variations seen in RX J0720.4— 3125 
during its spin period. However, to explain the change in the 
phase lags over time seen in RX J0720.4— 3125, requires a 
slow apparent movement of the spots relative to each other 
as seen from the distant observer. 

One possible explanation would be free precession of 
the neutron star as suggested by de Vries et al. ()2004p 
and put forward by H06. In this case the viewing geom- 
etry changes and no physical movement of the hot spots 
on the stellar surface is required. This model is discussed 
in detail in H06, but as already discussed above, can only 
be verified when the periodic behaviour is confirmed. The 
simple, free precession model discussed in II06 was able 
to account (at least qualitatively) for the long-term evolu- 
tion of the spectral parameters over the period examined 
(May 2000-November 2005) and also by invoking two non- 
antipodal caps for the peculiar intensity-hardness anti cor- 
relation with spin phase. Although this was not discussed 
in II06, the model predicts a cyclic variation in the phase 
residuals due to the change in the pulse shape (of the phase 
of the light curve maximun in particular) over the preces- 
sion period. Since the spin-phase dependence of the hard- 
ness ratio changes along the precession cycle (this can be 
also seen in the lower right panel of Fig. 5 in H06), one 
expects that the phase residuals and their variation over 
time are not the same in different energy bands. Indeed, 
as we checked, the free precession model outlined in II06 is 
able to reproduce (again, qualitatively) the properties of the 
observed residuals (see Fig. [2|. The highly non-sinusoidal 
pattern of the spectral parameters, which clearly emerges 
from the latest data (see Fig. 



Ill, however, may be, an 



issue. In order to explain the fast change around MJD = 
52800 days one needs to invoke the coming into view of a 



rather small "spot" with angular size '-^ At/Pp- 



10° 



10 



Hohle et al.: Spectral and temporal variations of RX J0720.4— 3125 





5.3 
MJD [days] 



i.15 5.2 



5.25 



5.3 
MJD [days] 



5.35 



5.4 



5.45 
X 10* 



Figure 11: Long term variations of the pulse phase averaged 
spectra obtained from the eleven EPIC-pn observations in 
full frame mode with a thin filter (dots) . Open circles show 
data obtained from the other EPIC-pn observations (see 
Table [ij. Error bars denote 2.7a confidence level. 



Figure 12: Phase resolved variations of equivalent width ver- 
sus temperature. Observations already used for Fig. 6 in 
H06 are marked in grey. The four new EPIC-pn data are 
coloured. The error bars denote 2.1a confidence level. 



where we take At ~ 500 d as the timescale over which 
the rapid change occurred and Pprec ~ lOyr. However, the 
presence of such a small emitting region produces a pulse 
shape which is not in agreement with the observed one. It is 
still under debate whether free precession can continiue for 
more than a few hundred cycles (Link ,2006.) . If this is the 
case and RX J0720.4— 3125 is precessing, then the preces- 
sion needs to be powered by some mechanism. Precession 
caused by an unseen sub-stellar companion can probably be 
excluded on the basis of the current upper limit on its mass 
obtained by Posselt et al. ((2008)) from a search for orbiting 
objects around RX J0720.4— 3125 (see also the discussion 
on PSRB1828-11 by Liu et al. [2006|) . 

Alternatively, the regions on the surface responsible for 
the hard and soft emission may physically change their loca- 
tions or their emission properties. A slow rearrangement of 
the magnetic field may change the local angle between the 
surface normal and the field, influencing the emission pat- 
tern, and the heat transport in the crust. In such a picture 
the rearrangements do not need to be periodic, but a differ- 
ent physical mechanism is required to explain it. A glitch 
(vK07) may be one possibility. The occurrence of such a 
sudden event might be supported by the relatively fast 
change observed in th e sp ectral parameters around MJD 
— 52800 days (see Fig. 1 1 1 and the more gradual evolution 
afterwards which resembles a relaxation process. However, 
the EPIC-pn small window mode data from 2003 (satellite 
revolution 711) and in particular the RGS data indicate a 
more gradual temperature increase before MJD = 52800 
days (de Vries et al.. I2004p . The small window mode data 
(thick and medium filter) was not used here at all and only 
partly by HOG, because of a normalisation discrepancy be- 
tween small window and full frame mode in the EPIC-pn 
calibration. This affects the inferred black body radius, but 
not the temperature. It should also be noted that the basic 
characteristic behaviour in the spectral variations was not 
altered around this date. The changes in kT, radius and 
EW seen on long term time scales are also seen during the 
neutron star rotation, before and after the putative event 
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(Fig. 12 and compare to Fig. 4 in HOG) and the total X- 
ray flux stays nearly constant within a few per cent (H06). 
The amplitude of the EW variation during the spin period 
did not change over the years, only the mean level evolves. 
In the case of the temperature, both amplitude and mean 
level evolve on long time scale, but this evolution is visible 
in the EPIC-pn spectra during all observations, first slowly 
and then fast with a large change around MJD = 52800 
days. Therefore, we regard a sudden event like a glitch as 
an unlikely cause for the (more gradual) spectral variations. 

In conclusion, the timing and spectral properties of 
RX J0720.4— 3125 are both consistent with the existence 
of at least two hot spots with different temperature on the 
surface of the neutron star. Whether the long term changes 
are periodic requires further monitoring. A simultaneous 
modeling of spectral and timing parameters is in progress 
and should eventually yield a more detailed temperature 
map of this unique isolated neutron star. 
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